THE SOUND OF AN EVOLVING FLOATING SCULPTURE 



BENJAMIN SEIBOLD AND YOSSI FARJOUN 

Abstract. Commissioned by MIT's in-house artist Jane Philbrick, we evolve an abstract 2D 
surface (resembling Marta Pan's 1961 "Sculpture Flottante I") under mean curvature, all the 
while calculating the eigenmodes and eigenvalues of the Laplace-Beltrami operator on the re- 
sulting shapes. These are then synthesized into a sound-wave embodying the "swan song" of the 
surfaces as the evolve to points and vanish. The surface is approximated by a triangulation, and 
we present a robust approach to approximate the normal directions and the mean curvature. 
The resulting video and sound-track were parts in the Jane Philbrick's exhibition "Everything 
Trembles" in Lund, Sweden, 2009. 



1. Introduction 

This paper presents a project in the areas of Computation and Apphed Mathematics. With 
both authors being researchers in Apphed Mathematics and Scientific Computing, this does not 
sound unusual. However, the area of apphcation is quite unusual (at least for the authors), namely: 
visual art. Employing and devising a variety of methods in computational geometry, numerical 
linear algebra, and spectral methods, two pieces of art — an animation and a sounds piece — were 
created. They were part of the exhibition "Everything Trembles" that ran at the Skissernas 
Museum (Lund, Sweden) in fall 2009. In this paper, we present the mathematical modeling that 
was required to formulate the artistic problem, the computational approaches undertaken, and 
the steps from computational results to creating pieces of art. 

In March 2009, through our colleague Enno Lenzmann (now University of Copenhagen), the 
visual artist Jane Philbrick (MIT Center for Advanced Visual Studies) approached the authors, 
with the following idea. She was preparing an exhibition fully devoted to Marta Pan's "Sculpture 
Flottante I", which is an exhibit in the Kroller-Miiller Museum (Otterlo, Netherlands). Her goal 
was to bring the sculpture to "life" in drawings, painting, and sculptures, as well as in modern 
visual art techniques, such as photographs, video, and sound installations. The specific question 
that she addressed us with was whether one can predict the sculpture's "mathematical future", 
by evolving and deforming it, as well as its characteristic sound. 

Not unexpectedly, our first step was to "model" the artist's idea into more precise mathematical 
concepts. Our answer, then discussed and agreed on by the artist, was to let the sculpture (as 
defined by its 2D surface) evolve under a mean curvature flow (thus answering the question about 
the future of the sculpture), and to create sound-form based on eigenfunctions of the Laplace- 
Beltrami operator on the evolving surface. Both ideas move the sculpture away from its rigid shape, 
and instead interpret it as a flexible membrane with tension, very much like a punctured soap- 
bubble-like surface. The mean curvature flow pulls the surface together, creates pinch-off events 
that break the object apart, and eventually leads to its complete vanishing. The Laplace-Beltrami 
operator acts analogous to the classical Laplacian, but on a curved surface. Its eigenfunctions are 
used to evolve wave forms in time that are in turn used to create a characteristic sound pattern 
that changes as the surface evolves under the mean curvature fiow. The two pieces of art created 
out of these mathematical concepts are: a video showing the evolution of the surface under mean 
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curvature flow, and a sound flle created out of the characteristic sounds of the surface as it evolves 
under the flow. 

This paper is organized as fohows. In Sect. [2j an outhne of the mathematical background 
relevant to this project is given. The "transformation" of the original sculpture into data that 
can be used for computational approaches is described in Sect. [3j The methods used and devised 
to compute the mean curvature flow are presented in Sect.|4j and the approaches to calculate the 
characteristic sounds of the surface are given in Sect. [5] Finally, in Sect. |6j we describe how the 
computational results are assembled into pieces of art. 

2. Mathematical Background 

2.1. Curvature. In the following, we consider smooth, simple curves in and smooth, ori- 
entable, and closed surfaces in R^. At any point on the curve/surface, h is the outward normal 
vector. For a smooth curve, C, embedded in R^, the curvature at a point P has an absolute 
value equal to the reciprocal of the radius of the osculating circle, i.e. the circle that passes through 
P and a pair of additional points on the curve that lie inflnitesimally close to P [26]. Moreover, 
we use the convention that hz is negative if n points towards the center of the osculating circle, 
and positive otherwise. With this convention, a circle of radius r has a constant curvature ^ = ^. 

For a smooth surface, embedded in R^, multiple notions of curvature exist. At any point, P, 
on 5, one can consider all curves on S that pass through P, and flnd the corresponding curvature 
for each of them (using the concepts established above). Of those curvatures, the maximal, /^i, 
and the minimal, hZ2^ are called the "principal curvatures" of S at P. Their product, hZitZ2j is 
called "Gaussian curvature", and their average, = ^(^i + ^^2), is called "mean curvature" [26] , 

While Gaussian curvature as an intrinsic quantity is independent of the surface's (isometric) 
embedding, mean curvature is extrinsic, i.e. it depends on the embedding of the surface in its 
ambient space R^. Mean curvature is closely related to the physical effect of surface tension, 
which acts normal to the interface between two immiscible fluids and is proportional in magnitude 
to its mean curvature. The mathematical relevance of the mean curvature was flrst formulated 
in the Young-Laplace equation [27l [TT] , which describes the scenario in which surface tension is 
balanced by a capillary pressure difference sustained across the interface between two static fluids. 

When there is no pressure counter-acting the surface tension on an interface between fluids (or 
a thin fllm) will move it as to decrease to surface area. A mathematical model for such movements 
due to surface tension is mean curvature flow [6j , which is a geometric nonlinear partial differential 
equation [6 . According to this flow, a surface embedded in R^, or a curve embedded in R^, moves 
parallel to its normal vector with a velocity proportional to its curvature. Speciflcally, a point x 
on the surface/curve moves according to the law of motion 

X — —tz{x)n{x) , (1) 

where is the (mean) curvature, and h is the outward pointing normal vector, both evaluated 
at X. The negative sign is due to our convention that curvature is positive where the surface is 
convex, and negative where it is concave. 

Curvature flow of curves in R^ exhibits a relatively regular behavior. Any simple, closed curve 
shrinks continuously, and at a steady rate, until it eventually vanishes into a single point 0. Due to 
the Gauss-Bonnet theorem [8], the total curvature over the curve is an invariant, i.e. hvdS = 27r, 
and thus the rate of change of area enclosed by the curve equals j^A{t) = x-hdS = — hidS = 

— 27r. Hence, the time of the curve's vanishing is T = where A{0) is the area enclosed by the 
initial curve. To the authors' knowledge, it is yet an unsolved problem to express the exact point 
of the curve's vanishing in terms of geometric properties of its initial conflguration [23 . 

Mean curvature flow of surfaces in R^ is far less regular than curvature flow of curves. First, the 
Gauss-Bonnet theorem does not apply to mean curvature (it applies to Gaussian curvature, which 
is not considered here). Second, a surface that moves under mean curvature flow can separate 
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into multiple disjoint components. Separation happens in a "pinch-off" event, in which the neck 
of a dumbbell-shaped part of the surface shrinks to radius zero. Since mean curvature flow is 
closely related to surface tension of physical fluids, this pinch-off phenomenon is well-known in 
fluid dynamics as a Plateau-Rayleigh instability, in which a falling stream of fluid breaks up into 
smaller droplets [l8l[T9]. In fact, in the evolution of "Sculpture Flottante I", we observe such 
pinch-off events. 

2.2. Laplace-Beltrami Operator. Since we are also interested in the "characteristic sound of 
the sculpture", we need to consider modes of vibration of a surface in R^. In Sect, [sj we use the 
surface wave equation 



as a model for sound wave propagation on the surface, where L is the Laplace-Beltrami operator, 
which is the generalization of the standard Laplace operator to curved surfaces. It is defined as 



where g is the metric tensor associated with the surface, and \g\ = \ det{gij)\ is the absolute value 
of its determinant. Since the Laplace-Beltrami operator encodes information about the surface's 
metric tensor, there is a close relation between it and mean curvature, namely (with our sign 
convention for curvature) 



where x is the position vector of a point on the surface, and H = tzh is the mean curvature normal 
vector [25 at the same position. As a consequence, the equation of mean curvature flow ([T]) can 
also be written in terms of the Laplace-Beltrami operator as 



which illustrates the similarity of mean curvature flow to a diffusion process: ripples in the surface 
shape are smeared out at a rate that is faster for ripples with higher spacial frequency. 

The surface wave equation ([2| is a linear partial differential equation, whose solutions can be 
conveniently expressed in terms of eigenfunctions of the Laplace-Beltrami operator. Clearly, these 
eigenfunctions relate to important questions about the vibration of objects, such as Mark Kac's 
fundamental question "Can one hear the shape of a drum?"[^ Moreover, the spectrum of the 
Laplace-Beltrami operator plays an essential role in many applications that go beyond the study 
of sound waves. Examples include signal processing and mesh design [24], shape matching [21], 
geometry processing [13j, shape analysis and image segmentation |2Q] . 

2.3. Computational Geometry. In order to implement the mathematical concepts described 
above, specifically mean curvature flow and eigenmodes of the Laplace-Beltrami operator, on a 
computer, these geometric operators must be approximated on a surface that can be described 
by a finite amount of data. We consider the surface be given by a regular triangulation. Hence, 
we must employ approaches to approximate surface normals, curvature, and the Laplace-Beltrami 
operator on such a discrete geometry. The first approaches of approximating differential operators 
on triangulations go back to Delaunay [4 . Since then, a wide variety of methods has been pre- 
sented. However, to the authors' knowledge, the quote "Despite extensive use of triangle meshes 
in Computer Graphics, there is no consensus on the most appropriate way to estimate simple 
geometric attributes such as normal vectors and curvatures on discrete surfaces" [15 still holds 
true nowadays. One key challenge is to devise efficient approximation approaches that work in a 
robust fashion under general circumstances. Below, we outline a few representative approaches 
for the approximations of curvature and the Laplace-Beltrami operator. 



utt = jCu 



(2) 




jCx = —2k{x) , 



(3) 



X — 2 «^ 7 
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Many approaches to approximate the Laplace-Beltrami operator at a vertex of a triangulation 
are based on so-cahed "cotan formulas" . Let V denote the set of vertices, and E the set of pairs of 
vertices that are connected by an edge. Assume that a function u is given on the vertices. Then 
>Cii at a vertex Xi is approximated by 

Cu{xi) ^ ^ Wij{u{xj) -u{xi)) , (4) 

where Wij are appropriate weights. A common choice for the weights [2 is given by 

Wij = \{coiaij ^coi a ji) , (5) 

where aij and aji are the angles that are opposite to the edge (xi^xj)^ in the two triangles 
adjacent to this edge (see, for example. Fig. 1 in [2]). Other, similar formulas, albeit with different 
normalization constants, have been studied in the literature [5, 15 . 

Due to relation ([3|, an approximation rule for the Laplace-Beltrami operator induces an ap- 
proximation to the mean curvature normal vector ^ = /^n, by applying the formula to the surface 
coordinate vector x. For instance, formula Q leads to 

^ Wij{xi — Xj). (6) 

XjeV:{xi,Xj)eE 

In words: an approximate curvature normal vector at a vertex Xi is obtained as a weighted sum 
of all adjacent edge vectors that point into Xi. In many cases, cot an formulas yield good quality 
approximations to mean curvature, but not always: for example, this method does not perform 
well on extremely pointy configurations [5, Sect. 4.1]. Another shortcoming of cotan formulas is 
that they do not reveal the full curvature tensor. Examples of alternative approaches include the 
"Functional Laplacian" [l], which is based on Green's functions, and approximations using the 
Willmore energy of surfaces [9]. An alternative way to approximate curvature is least-squares 
fitting for principal directions [15]. 

In this project, we employ the classical formula Q for the computation of the spectrum of 
Laplace-Beltrami operator (see Sect.|5|, and a non-classical way to approximate curvature (based 
on least-square minimization) to achieve a stable and robust approximation of the mean curvature 
fiow. This latter approach is described in detail in Sect. |4.1[ 

It should be remarked that there are other computational methods for mean curvature fiows 
that are not based on explicit representations (e.g. triangulations) of the surface. An important 
example is the level set method [16 , which represents the surface implicitly as the zero contour of 
a level set function cj). One then can formulate a partial differential equation (in the ambient space) 
for so that its zero contour moves according to mean curvature fiow [23l[3]. Advantages of the 
level set approach are the conceptual simplicity of solving a PDE in the ambient space, and the 
automatic treatment of topology changes. Disadvantages are the extra problem dimension (due 
to moving to the surface's ambient space) and the loss of sharp surface features due to numerical 
viscosity. In addition, if the original configuration is given as a triangulation (as it is in our case), 
level set approaches would require back-and-forth transformations between explicit and implicit 
representations, each of which would incur errors. 

3. Transformation of the Sculpture into a Computer Model 

The "Sculpture Flottante I" is a permanent exhibit in the sculpture garden of the Kroller-Miiller 
Museum (Otterlo, The Netherlands). It was created in 1961 by Marta Pan (1923-2008), a French 
sculptor of Hungarian birth. The sculpture, shown in its installation view in Fig. [T] is made of 
fiberglass reinforced polyester. It consists of a base which is fioating on a pond, and a top which 
is shaped like a cut-open umbrella and is connected to the base by a thin shaft. 

In order to admit computational tasks with the sculpture's surface (as well as its visualization), 
a computer model needed to be generated. Since a direct scanning of the original sculpture turned 
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Figure 1. Marta Pan's "Sculpture Plot- Figure 2. Creation of Figure 3. Laser scanning of the 
tante I" at the KroUer-Miiller Museum|^ a clay model[^ model|^ 

out not possible, a clay replica of the sculpture was created (shown in Fig. [2|, based on a series 
of photographs of the "Sculpture Flottante I" . The clay model was then laser scanned (shown in 
Fig.[3| resulting in a stereolithography CAD file (.stl format), from which a triangulation of the 
2D surface that resembled the sculpture could be extracted. This triangulation, which consists 
of more than 500,000 elements, was then transformed into a format that allowed us to perform 
computational geometry operations, as described in Sect. |2.3l Sect.|4j and Sect. [5] 

As can be seen in Fig. |2j the true sculpture's top and bottom are connected through a very 
thin shaft. Mean curvature flow contracts this shaft almost immediately, resulting is a pinch-off 
that separates the top from the bottom. We choose the time of this pinch-off as the starting time 
of the animation, described in Sect. [6j 



4. Evolution of the Geometry under Mean Curvature Flow 

The computational task is to let the initial 2D surface evolve under mean curvature flow ([T]), 
until it has vanished completely. The basic idea is to advance equation ([T]) forward in time using 
an explicit forward Euler scheme. This amounts to the following update scheme. In each step: 

i. Given the triangulation, evaluate the mean curvature normal vector in the right hand of ([T]). 

ii. Move each vertex of the triangulation by At times this right hand side vector, where At is an 
appropriately chosen time step (see below). 

iii. If required, modify the triangulation (see below). 

In order to make this basic approach work, several technical challenges must be addressed. 

First, the mean curvature normal vector on the vertices of a triangulation must be approximated 
in a robust fashion, such that motion under ([T]) results in a stable scheme. It must be stressed 
that the approximation of the curvature normal vector for the sake of an evolution equation is 
significantly more demanding than a single evaluation of these differential quantities on given 
data. As an analogy, consider the approximation of the first derivative of a function using central 
differencing. When done once, this is a good idea. However, when used with forward Euler time 
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stepping to approximate the linear advection equation, an unstable numerical scheme results. The 



methodology employed here is described in Sect. 4.1 



Second, as particles move, the topology of the surface may change (in a pinch-off event). The 
computational approach must be able to mimic such topology changes. This is described in 
Sect. [O 



Third, as the surface shrinks, particles tend to move closer towards their neighbors, and as they 
get closer the largest admissible time step shrinks. To guarantee that a sufficiently large time step 
can be taken, we merge neighboring vertices if they become too close, and then check for possible 



resulting topology changes. This aspects is described in Sect. |4.3 

Finally, the aspect of computational efficiency must be addressed. The maximum admissible 
time step is governed by the smallest features of the complete geometry. As a consequence, in 
most other parts of the surface, this time step is unnecessarily small, and thus a full update of 
the curvature normal vector everywhere would be prohibitively costly. Our remedy, described in 
Sect. |4.2[ is a mechanism by which the normal vector and curvature are only calculated where 
and when they truly need to be re-calculated. As long as they do not require an update, the 
last calculated values are used. Effectively, this results in a scheme with spatially and temporally 
adaptive local time steps. Since the computation of the differential quantities is a costly part of 
the computation, this added adaptivity reduces the computation time dramatically. 



4.1. Approximation of Normal Vectors and Mean Curvature. As described in Sect. 2.3 
there is a wide variety of methods to approximate surface normal vectors and curvature, some 
of the most popular approaches being the "cotan formulas" Approximations based on these 
formulas can be used to compute mean curvature flows, at least for short times, as it is sufficient 
for surface smoothing (see [5 for an example). For the task of evolving a surface all the way to its 
vanishing, and including topology changes, we have found that an alternative approach, described 



below, is more robust. This approach finds the normal vectors in a first step (see Sect. 4.1.1), and 



then finds the curvature tensor using a least-squares fitting approach (see Sect. 4.1.2). While this 
two-step approach turns out to yield a robust mean curvature flow evolution, we are admittedly 
overdoing things a bit in having access to the normal vector and the full curvature tensor. As a 
consequence, the approximations presented here have the potential to be used for geometric flows 
that are more complicated than mean-curvature flow — a topic that we do not follow up with in 
this paper. 

4.1.1. Approximation of normal vectors. To approximate the normal vector at a vertex vq, we look 
at two types of averages and chose the larger between them. In the following the center vertex vq 
is assumed to be at the origin as this simplifies the notation. Moreover, we consider the i direct 
neighbor vertices vi^ . . . be ordered counterclockwise around the center vertex when the shape 
is viewed from the outside. This is possible for all orient able surfaces. The choice of vq as the 
origin implies that the vertices Vj are also the edges connecting the vq to the vertices. 

• The average face normals of adjacent faces around the vertex, weighted by the angle at 
the vertex (see Fig. |4|: 

E} = l^i^(0,i,j + 1) 



Face normal average 
where 



arcsm 



is the angle between the edges j and j + 1, and 

_ Vj X Vj^i 

is the normal to the face that is between the j and j + 1 edges. 
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Figure 4. The face- aver aged normal vector. 



Figure 5. The edge- averaged normal vector. 



• The average of normalized adjacent edges^ pointing toward the center vertex, — Vj, weighted 
by the sum of the two angles at the vertex and that edge (see Fig. [5|: 



Edge direction average 



(7) 



Note that ^ is similar in flavor to the cotan approximation ([6| to the curvature normal 
vector. 

Of these two averages we use the one that is larger in magnitude, since this would imply that 
the parts being averaged agree on the normal to a greater extent. The reason we are using two 
different methods is that we need to flnd an accurate normal direction in both extremal cases: 
(a) when the faces are part of a large, nearly flat domain, and (b) when neighboring triangles are 
at the knife's edge of a surface, or are part of a very narrow spike. In the friendly case (a), the 
face normal average yields a good normal direction, since the face normals are almost parallel. In 
contrast, the edge direction average possesses only a small component in the direction normal to 
the surface, and is henceforth not used. In the case (b), when the center vertex is on a tip of a very 
sharp, needle-like, formation, the face normal average is very small (and thus discarded), while 
the edge direction average is nearly a unit vector in the same direction as the sharp formation. 

When evolving a surface under mean curvature flow, in the vast majority of times the close- 
to-flat case (a) is observed, as the regularizing flow tends to smooth the surface. However, the 
proper treatment of the "spiky" case (b) is still required to obtain a robust treatment of the initial 
conditions, and of pinch-off events. Whichever type of average is used, the resulting vector is, of 
course, normalized before being used as the surface normal. 



4.1.2. Approximation of the curvature tensor. Having found a normal vector h at the vertex vq^ 
we now need to approximate the curvature. To do so, we introduce a local coordinate system, 
where z = h is that approximate normal, and x and y span the plane n^. The following formulae 
are written with regards to that coordinate system, where vq is the origin and the points (x, z) 
refers to vo-\-xx-\-yy-\- zz. In these coordinates the neighboring vertex vj is located at (xj, yj^ zj). 

We now approximate the curvature by fltting a quadratic polynomial 

a -\- x^ ^ /) ^ 

to the triangulated surface IJ}=i('^o = where, for convenience of notation, we deflne 

Vi-^i = vi^ and furthermore x = (x^y)^. To approximate the shape constant (i.e. the curvature 
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tensor and the constant shift) of the surface, we look to minimize the functional 



/triangle j 



d e 
e f 



z{x,y) dxdy , 



(9) 



where w = (a, e, f)^ is the vector of shape constants to the surface at a given vertex, and 
z{x^y) is the z- value of the triangulation at the point {y). We re- write ^ using a vector of 
functions, Sf^x^^ = (l, x^, 2xy^ y'^) , as 



^ J J triangle j 



s^x w — z{x^ y) ) dx dy 



In each triangle j we use barycentric coordinates a, /3 and 7 (using the three vertices of the j— th 
triangle) where a + /3 + 7 = 1, to identify a point x inside the triangle: 



Xj Xj-^i 



Due to the choice ^0 = 0, the value 7 does not appear in the expressions. Changing the integration 
variables to a and /3, we obtain 



i=i 



Xj 




nl nl — a / 




%+i 


^0 ^0 \ 



dp da , 



where the determinant of the Jacobian enters due to the change of variables. 

From the definition of s one can verify that applying a transformation Xj = yjXl ) to the 
argument of s is equivalent to right-multiplying by a matrix: 



/I \ 

x^j 2xjyj y] 

XjXj^i Xj^iyj ^ Xjyj^i yjyj+i 



(10) 



The non-trivial lower 3x3 block in the matrix in (10) is 5'^(Xj), the second symmetric power of 
Xj. Letting a = (^), and Zj = {z%-i ), the minimization functional can be written as: 



^[w]=^|X,| 



1 

^0 



w — a Zj 



dj3da . 



To integrate, we expand the power as a product of a vector with its transpose: 



EM=^\X, 

3 = 1 



^0 



w S - S(^s-SjW — 2w S • Sc^a Zj + Zj aa Zj dp da. (11) 



In (11) the only terms that depend on a are s^s^, s^a^ , and aa^ . These terms' integrals are 



T 



1 

180 



/90 15 15 15\ 

15 6 3 1 

15 3 4 3 

\15 1 3 6 / 



-"T _ 1 



/lO 10\ 

3 1 

2 2 

VI 3/ 



_ X l"^ 1 

aa — 24 I 2 2 



which we label B^C respectively. 
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To find where the minimum of E[w] is attained, we replace the integrals in (11) by A, 5, and 



C and differentiate with respect to w. This gives rise to the linear system for the minimizing w: 

^ |X, I EjAE, I w = ^ |X, I EjBzj . (12) 

The mean curvature is now given by z^: = = (0,^,0,^) -w. Note that according to our 
approximation, the surface does not need to pass through any of the vertices vq^. . . ^Vi^ but the 
tangent plane at vq is perpendicular to n. It is in principle possible to allow the fitted form to 
have a non-trivial tangent plane by including a term 6x + c?/ in ([8|, however, we have not found 
this improves the quality of the results. 




4.1.3. On the non-robustness of vertex-based minimization. The previous section might seem ex- 
cessive effort for fitting a quadratic form. A more straightforward approach might be to simply 
minimize the distance of the form to the vertices: 

It has been observed (for example by Meyer et al. [15 ) that this minimization can give poor 
results. Mathematically, minimizing this functional results in a linear system that may be under- 
resolved. For example, if vq, • • • , ^4 comprise a rectangular pyramid, there are two trivial quadratic 
forms (and therefore an infinite number of them) that attain E[w] = but have different mean 
curvatures. A vertex with three neighbors will always have a form that attains E[w] = 0, and 
the resulting quadratic forms often do not resemble the surface they are intended to approximate. 
This effect is drastic, and can result in a mean curvature of the opposite sign as the correct one, 
unstably causing a spike to emanate from the surface and grow indefinitely. 

While we have no proof of robustness of our approach based on £^[w], during the evolution 
of the "Sculpture Flottante I", with its 500,000 vertices to start, through merges and topology 



changes, until its demise not even once was linear system (12) poorly conditioned. By minimizing 
over the triangles (and not only the vertices) we also provide a robustness with regards to refining: 
if the triangulation is refined by adding a vertex on the edge which connects two neighbors of vq^ 
our functional E[w] remains unchanged while E[w] changes. 



4.2. Staleness. When a triangulated surface evolves under mean curvature, the curvature at some 
vertices changes very quickly, while at other it changes slowly. Calculating the curvature at all the 
vertices at small time intervals can be computationally expensive, however, to get stable results 
we must calculate it at short time intervals for vertices where curvature changes quickly. 

To get around calculating the curvature for all the vertices at the rate set by the vertices who 
need the most frequent update, we recalculate the curvature for some vertices more frequently 
than for others. This is done by estimating for how long the curvature calculation remains relevant 
("fresh") at the same time that it is calculated. The "expiration time" we use for a vertex vq is 
given by 

Atexp = C, min 

or, in words, some constant (we used Cs = 0.5) times the length of the minimal adjacent edge, 
divided by the velocity at which the the vertex moves, k,. 

To avoid taking exceedingly small steps (stepping from one expired vertex to the next), once 
the curvature needs to be calculated at a vertex, we also recalculate all the vertices whose data is 
about to expire in the next time step. Then the largest step that can be taken (based on staleness) 
is performed. Of course, we also limit the step size by standard accuracy considerations, and make 
sure that we reach certain points in time to produce an output. 
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Figure 6. A topology change is initiated when two of the neighboring vertices vo, vi, and V2, are to 
be merged into a vertex vi, but the triangle (vo, vi, V2) (the shaded triangle), is not part of the surface. 
After the merge, the resulting surface is not homeomorphic to the 2-disk around vq and vi. To resolve 
this, vo and then vi are duplicated (generating Vq and Vi) and then connected to a subset of their original 
neighbors. 



4.3. Removal of Vertices. Since the size of the largest time step we are allowed to take is limited 
by the proximity of neighboring vertices, we take care to remove vertices that are too close to each 
other. We used three conditions as triggers to remove vertices: 

• small edge, or 

• small angle, or 

• small face. 

To remove a vertex, we merge two close vertices into a new one whose location is an average of 
the two. The "freshness" of the neighboring vertices (and the new vertex) is set to zero so that 
the curvature there be re-calculated immediately. 

Special care needs to be taken to identify a required change in topology, and then make it 
happen. A change in topology is required when the neighborhood of a vertex is not homeomorphic 
to the 2-disk. As visualized in the stylized Fig. [6j this can happen when two neighbors, vi and 
of a vertex, vq^ are merged, to vi^ but the surface does not contain the triangle {vq^ vi^ V2). This 
happens when the triangle forms the neck of a part that is about to pinch off. To identify this 
event, we examine the set of '^o's neighbors, ['^i, . . . , Vi] and attempt to order it so that for all j, 
Vj and Vj-^i are neighbors. If this is not possible (because some vertex is needed in more then one 
place), there is a topology change. In this case, we find a subset of the vertices for which this is 
possible, and to each (perhaps not disjoint) set, we "connect" a unique copy of vq. After splitting 
Vq, the same operation will be needed for vi, as shown in Fig. |6] These types of topology change 
treatments are similar to the ones presented by Mayer [14 , with the difference that in said paper 
surface diffusion flow is considered, and not mean curvature flow. It should also be remarked 
that unlike other types of geometric flows, for mean curvature flow the insertion of vertices is not 
required. 
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5. The Characteristic Sounds of the Surface 



As described in Sect. 2.2, we use the surface wave equation ([2| to characterize the "sound of 



the sculpture" as it evolves under mean curvature flow. Since we are dealing with a surface that is 
evolving in time, S{t)^ we also have a time-dependent family of Laplace-Beltrami operators. We 
denote this dependence on the surface by a subscript, i.e. the function u{t^x) which we use to 
characterize the "sound of the sculpture" , as it evolves under mean curvature flow, satisfies 

Utt(t,x) F^Cg^^^u{t,x) = for ah t > 0, and x G S{t). 

To simplify, we assume that the sound-generating vibrations have a much shorter time scale than 
that of the evolution of the surface itself, i.e. ^ 1, and we allow ourselves to use a quasi-static 
approximation by which the surface itself is assumed to be constant between time steps. At every 
time step, we project the solution onto the evolved surface and continue with the time evolution. 

For a fixed surface, = S{t^)^ solutions to 

Uu{t,x)^F^C^^u{t,x) = (13) 

can be found using separation of variables: 

u{t,x)=e^'^^^'^{x) , 

where 

-C^^(I) = X(I). (14) 
Note that in the following, we consider the eigenvalues of the negative Laplace-Beltrami operator. 



since this is a positive semi-definite operator. The general solution to (13) can be written using a 
linear combination of such simple terms: 



oo 



fe = l 

where ((/)^, A^)^-,^ are all the eigenfunction/ value pairs of C^^- 

After evolving on the fixed surface S'^ for a time At, we project it and its time derivative, 
u^^ onto the surface S^^^ . The projected functions are used as initial condition for the evolution 
for another time step of At. From the resulting coefficients ^ and the eigenvalues we create 
the characteristic sound-wave of the evolving surface. 

5.1. Eigenfunctions and Eigenvalues. To approximate the Laplace-Beltrami operator C^ri on 
the surface S{t^), we employ a "cot an formula" Q with the specific weights (jsj, as it is given 



in [2]. This discretization replaces the infinite dimensional eigenvalue problem (14) by a matrix 
eigenvalue problem 

where (p = {(p{vi)^ . . . , (p{v\s^\))^ ^ and the entries of the matrix are given by the weights in the 
cotan formula, as 

a'^i = ^ Wij \/i and a^j = —Wij yi^j, (16) 

3 

which, of course, depend on the specific triangulation of at time . Consequently, the matrix 
is a large (initially about 500, 000 x 500, 000), sparse, symmetric matrix, and we need to find its 
smallest A/" = 50 eigenvalues and corresponding eigenvectors. We achieve this task by an implicitly 
restarted Arnoldi method implemented in Matlab's eigs. The symmetry of A arises from the 
self-adjointness of C and from the symmetry of the cotan formula. It implies that the resulting 
eigenfunctions are orthogonal (under summation over the vertices). They are normalized so that 

In this notation, (p^ denotes the k-th eigenfunction of the discrete Laplace-Beltrami operator. 



discretized according to ( 16 ) on the triangulation at time t = f^' 




The temporal evolution of the N lowest eigenvalues of the sculpture's top part, as it evolves 
under mean curvature flow, is shown in Fig. [7[ One can see how the frequencies increase as the 
surface shrinks. Moreover, one can observe the discontinuity in the eigenvalues due to the pinch-off 
event (shown in Fig. |9j between panels (c) and (d)). The eigenvalues that equal zero correspond 
directly to the connected components of the surface. The brief appearance of a second connected 
component after the pinch-off is echoed by a short segment in the "zero eigenvalues" zone. 



5.2. Initial Condition and the Evolution During a Time Step. Equation (13) requires 
initial conditions for both u and Ut. For all n > we use the projection from the end of the 
previous step as the initial conditions for the current step. The initial condition at n = is chosen 
to attain 1 on a single vertex (chosen manually), and to vanish on the rest. The time derivative, Ut, 
is chosen to vanish identically at n = so that we can visualize this initial condition as "plucking" 
the mesh at one vertex. 

At the beginning of each step the values of ^ are found from the values of u'^ and u'^ by 
solving a sequence of small linear systems: 




This system is found from ( |T5| ) by differentiating with respect to t and then taking the inner 
product with the eigenfunction (/^J^. 

Once the coefficients of the solution at time t = have been found, the solution and its time 



derivative are evolved to t = t^^^ using the general solution (15): 

N 



N 
k=l 

where we write u and not u to signify that it was not yet projected on S{t'^~^^). 
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5.3. Projection. Before we can take another time step we must project the function and its time 
derivative from the surface S{t^) onto S{t'^~^^). During the evolution, we keep track of the merging 
of vertices, so that for every vertex v^~^^ at time we can identify a "parent vertex" vf" at 
time with the property that the evolution of the surface takes vf" onto v]^~^^ by evolution and 
merges. This "parent image" is used for the projection: the value of a function /(^^^^) is taken 
to be fivj"). 

5.4. The Wave- Form. Once we have found the coefficients ^ for the evolution of the sculp- 
ture until it vanishes, we construct the characteristic sounds of the evolution. We approximate 
continuous functions Ak{t) and Ck{t) by a linear interpolation of the and the ^. Using these, 
we generate the waveform W{t) as the sum 

k=l ^ ^ 

which is then sampled at a high resolution to create a .wav file. 



6. The Transformation of the Computational Results into Pieces of Art 

Having created the computational results from the mean curvature flow of the sculpture's surface 
(see Sect.[4|, and having an approach to compute its characteristic sounds (see Sect.[5|, we now 
describe the crucial steps of assembling these data into pieces of art. 

6.1. The Video. A video is created that shows the evolution of the sculpture's surface under mean 
curvature flow, from its initial state, until it vanishes completely. The triangulation is rendered 
using a Phong shading illumination model [17j, involving ambient, diffuse, and specular reflections 
of two light sources. At any instance in time, each connected component receives a subtle coloring 
based on its three lowest (non-constant) eigenfunctions, i.e. the ones that correspond to the three 
smallest, positive eigenvalues (there is a zero eigenvalue for every connected component, corre- 
sponding to a constant eigenfunction on that component). Specifically, the lowest, second- lowest, 
and third-lowest frequency eigenfunctions affect the red, green, and blue of the surface's RGB 
channels, respectively. As described in [131 120], ^he resulting coloring of the surface's connected 
components adds some characteristic geometric information to the visualization. 

In discussion with the artist Jane Philbrick, a number of artistic features are added. First, 
the evolving surface is placed on top of a blue plane, and a mirrored, blue-shaded copy of the 
surface is placed below the plane, thus imitating reflection in water, as it would occur with the 
original "Sculpture Flottante I" (see Fig. [T]). This water-like effect is further amplified by adding 
small ripples on the reflection. Second, an animation sequence is added to precede the actual 
mean curvature flow, in which the digital version of the "Sculpture Flottante I" is introduced to 
the viewer. In this sequence, the sculpture is initially approached by a moving camera, and then 
rotated to allow a 360° view of its shape. Third, a counter is added that clocks the pseudo-time 
that has passed from the onset of the mean curvature flow, providing a reference for the various 
time scales that the flow possesses: the fast contraction of high-curvature regions on the one hand, 
and the slow shrinking of almost flat areas on the other hand. A frame of the video is shown 
in Fig. |8] At the instance shown, the mean curvature flow has just started, and created a small 
separation between the top and the bottom pieces. A sequence of eight representative snapshots is 
given in Fig.[9j The initial configuration is shown in panel (a). Between panels (b) and (c), a small 
pinch-off event happens (see Fig. [t] for the corresponding behavior of the eigenvalues) . Between 
panels (f) and (g), the top piece vanishes, and finally in panel (h) the base piece is approaching 
the last stages of its existence. The full video can be downloaded on the project website [22j. 
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Figure 9. The evolution of the sculpture's surface under mean curvature flow. From top left to bottom 
right: (a) the initial configuration, and the surface after (b) 100 , (c) 125 , (d) 130 , (e) 200 , (f) 400 , 
(g) 700 , (h) 1500 time units. Between (c) and (d) a small piece is pinched-off the top part. 

6.2. The Sculpture Song. The approach described in Sect. |5] yields the temporal evolution of 
an initial "plucking" of a single point in one connected component of the surface, as it evolves 
under mean curvature flow. When, where, and how frequently to excite the pieces of the surface 
to create an interesting sound piece is an artistic decision. As such, we created a number of sound 
files, arising from different initial conditions and provided them to the artist, who then mixed 
these pieces together to create the final "Sculpture Song" . The resulting sound file is available for 
downloaded on the project website [22j. 
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